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Abstract 

We present numerical simulations of fully nonlinear drift wave-zonal flow (DW-ZF) turbulence 
systems in a nonuniform magnetoplasma. In our model, the drift wave (DW) dynamics is pseudo- 
three-dimensional (pseudo-3D) and accounts for self-interactions among finite amplitude DWs and 
their coupling to the two-dimensional (2D) large amplitude zonal flows (ZFs). The dynamics of 
the 2D ZFs in the presence of the Reynolds stress of the pseudo-3D DWs is governed by the 
driven Euler equation. Numerical simulations of the fully nonlinear coupled DW-ZF equations 
reveal that shortscale DW turbulence leads to nonlinear saturated dipolar vortices, whereas the 
ZF sets in spontaneously and is dominated by a monopolar vortex structure. The ZFs are found to 
suppress the cross-field turbulent particle transport. The present results provide a better model for 
understanding the coexistence of short- and large-scale coherent structures, as well as associated 
subdued cross-field particle transport in magnetically confined fusion plasmas. 
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It is widely recognized that the presence of large scale sheared flows [also referred 

to as convective cells (CCs) or zonal flows (ZFs)] is detrimental to regulating the cross- 
field turbulent transport in magnetically confined fusion plasmas. The ZF is characterized by 
poloidally and toroidally symmetric structure with radial variation, and the relative zonal 
flow potential fluctuation (in comparison with T e /e, where T e is the electron temperature 
and e is the magnitude of the electron charge) is much smaller than the relative zonal flow 
density perturbation (in comparison with the equilibrium plasma number density no). The 
large scale Zonal jets also occur in various planetary atmosphere, where they are nonlinearly 
generated by the Rossby waves % 15], and influence the atmospheric wind circulation a, l3|. 

In magnetically confined fusion plasmas, there exist free energy sources in the form of 
density, temperature, and magnetic field inhomogeneities, which are responsible for exciting 
the low-frequency (in comparison with the ion gyrofrequency) , short scale (of the order 
of the ion gyroradius or the ion sound gyroradius) DW-like fluctuations {{j, 9, 10]. The 
linearly growing drift modes interact among themselves and attain large amplitudes in due 
course of time. The Reynolds stress of finite amplitude DWs, in turn, nonlinearly generate 



convective cells (CCs) and sheared flows/ZFs 



11 



12. 



and modulational instabilities 



16, 



171 ] , via three- wave decay 



12], respectively. There are recent review articles presenting 



the status of theoretical and simulation works 1171, as well as experimental observations 

n 

|18] concerning the dynamics of DW-ZF turbulence system. Specifically, some numerical 
simulations l_7j] lend support to the experimental observation that the DW turbulence and 
transport levels are reduced in the presence of the sheared flows/ZFs. 

Recently, Guo et al. [3] used the governing equations of Ref. 12] for the DW-CC 
turbulence system to investigate the radial spreading of the DW-ZF turbulence via soliton 
formation. However, the authors of Ref. [19| completely neglected self-interactions among 
drift waves and zonal flows, which are very important in the study of nonlinearly coupled 
finite amplitude drift and zonal flow disturbances in nonuniform magnetoplasmas. 

In this Letter, we present simulation results of fully nonlinear DW-ZF turbulence systems, 
which exhibit the coexistence of drift dipolar vortices and a radially symmetric monopolar 
zonal flow vortex. The effect of the latter on the cross-field turbulent transport is examined. 
Our investigation is based on the governing equations for the DW-ZF turbulence systems 
that incorporate the Hasegawa-Mima (HM) self-interaction nonlinearity [20] in the nonlinear 
dynamics of the DW which are nonlinearly exciting CCs/ZFs. Furthermore, we also account 
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for nonlinear interactions among the CCs/ZFs and obtain the driven Euler equation for 
the dynamics of finite amplitude CCs/ZFs. The generalization of the governing equations 
for fully nonlinear DW-ZF turbulence systems is rather essential for the investigation of 
the formation of coherent nonlinear structures that control the transport properties and 
confinement of tokamak plasmas. 

We consider a nonuniform magnetoplasma in an external magnetic field zBq, where z is 
the unit vector along the z— axis in a Cartesian coordinate system and Bq is the strength of 
the homogeneous magnetic field. The density gradient dn^/dx is along the In the 

presence of the finite amplitude low-frequency (in comparison with the ion gyrofrequency 
u ci = eBQ/rriiC, where m ; is the ion mass and c is the speed of light in vacuum) electrostatic 

DWs and ZFs, the perpendicular (to z component of the electron and ion fluid velocities 

13| are, respectively, 

« x V0 - -£—z x V(n e T e ) = u d EB + u d De , (1) 

-DO -DO^e 

u* ± » (c/B )z x = u| B) (2) 

uf ± « u d EB + u d Di - (Jt + u in - 0.3^p, 2 Vl + u| fl • V + i4 ■ V ± (3) 

>| B • V) V ± V + « s - V) V ± 



and 



B u c 

^ + z/ m - 0.3^Vi + (« B • V)Vj 



(4) 



where the superscripts d and z represents quantities associated with the DWs and ZFs, 
respectively, and ip are the electrostatic potentials of the DWs and ZFs, respectively, n e 
and iii are the electron and ion number densities, respectively, u Ei = [cj eB^n^i x V(njlj) 
is the ion diamagnetic drift velocity, Tj is the ion temperature, Vi n {vn is the ion-neutral 
(ion-ion) collision frequency, p, t = Vri/u} ci is the ion gyro-thermal radius, and Vn is the 
ion thermal speed, We stress that the self-interaction nonlinearities of the DWs and ZFs 
are retained in the fluid velocities (3) and (4), respectively. The angular brackets denote 
averaging over one period of the DWs. 



Assuming that (d/dt) + u EB ■ V 



<C u en , where v en is the electron- neutral collision fre- 



quency, we obtain from the parallel (to z component) of the electron momentum equation the 
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magnetic field-aligned electron fluid velocity u d z ~ (l/m e u en ) d {ecj) — TeTi^Jnoj /dz, where 
[n e — n ) <C n . We can now insert u d z into the electron continuity equation to obtain 







+ 



VI d 2 



dt Vpr, dz 2 



+ (u d EB + u z F , R ) • V 



EB 



n d el + u d EB ■ Vn + 



m P v en dz 2 



0. 



(5) 



where Vr e = {T e /m e ) l l 2 is the electron thermal speed and m e is the electron mass. Fur- 
thermore, substituting for the ion fluid velocity from (3) into the ion continuity equation we 
have 

~d 



Of 



+ [u EB + u 



l EB 



• V 



< + u% B ■ V(n Q + n\Z 



(6) 



n c 



B U C 



d_ 

dt 



^- + v in - 0.3u iiP 2 Vl + u d EB ■ V Vi0 + V ■ {u Di d ■ V) Vj 



n c 



B u c 



[u d EB ■ v)vi^ + {< B ■ v)v 2 



0. 



where the magnetic field-aligned ion dynamics has been ignored, thereby isolating the ion 
sound waves from our system. The ion density perturbation associated with the ZFs is 
n z a = (n c/B uj ci )V 2 L iJ. 



Equations (5) and (6), which govern the dynamics of collisional drift waves 2l| in the 



presence of zonal flows, are closed by assuming n d x 



n 



a 



rii, which is a valid approxima- 



tion in plasmas with S> uJ 2 i7 where uo pi is the ion plasma frequency. In the linear limit, 



without the ZFs, Eqs. (5) and (6) yield the DW frequency u>k = —k y c s p s / L n {l + k^p 2 ) and 
the growth rate 7fc(-C Uk), which are much larger than the damping rate v in + ^."ivak\p 2 . 
The growth rate is 7^ = v en {jj\k\j{jj 2 LH k 2 z {\ + k"j_p 2 ), where c s = (T e /rrii) l l 2 is the ion sound 
speed, p s = c s /uj ci is the sound gyroradius, cu LH = {uj ce uj ci ) l l 2 is the lower- hybrid reso- 
nance frequency, oj ce = eB /m e c is the electron gyrofrequency, L n = (dlnno/dx)^ 1 is the 
scale- length of the density gradient, and k = + zk z is the wave vector. 

The equation for the ZFs is obtained by inserting (2) and (4) into the electron and 
ion continuity equations, and inserting the resultant equations into the Poisson equation, 
obtaining the driven [by the DW Reynolds stress; the last term in the left-hand side of Eq. 
(7)] damped (by the ion-neutral collision and ion-gyroviscosity effects) ZF equation 



' d 



0.3u iiP 2 Vl + u z EB ■ V VlV 



u 



EB 



v)v: 



0. 



(7) 



For the collisionless DWs, we assume that \(d(j)/dt) + (u EB +u EB ) ■ V<j>\ <C (VrJ V en )Vj_0 
and z x Vn ■ V0 <C (u; ce /z/ en )n o |<9 2 0/<9z 2 |, and obtain from (5) the Boltzmann law for the 



electron number density perturbation n 2 el = noe0/T e . The latter can be inserted into (6) 
by assuming that n d iX = n^, so that we have fully nonlinear equation for the DWs in the 
presence of ZFs 

' d 



d(p _ c s p s d(f) 
dt L n dy 



Ps 



. + v in - 0.3z^Vl + — (1 + a) (z x V0) ■ V 

OT JDq 



+—(z x W) • V 

-DO 



2vv2 



o. 



where a = Ti/T e . 

We normalize the time and space variables by Ud-i and p s , as well as <ft an d by T e , 
and the collision frequencies by u C i. In the normalized units, we can rewrite (7) and (8) as, 
respectively, 



^ + ^ _ 0.3^aVl + (z x V^) ■ V 



VlV + ((zx V0-V)V^ 



0. 



and 



O0 _ 

ot L n dy 



d + ^ _ .3^xVi + (1 + <r)(z x V0) ■ V 

ot o; C j aj c j 

+(i x W) • V (0 - V 



V 2 



(9) 



(10) 



0. 



We have developed a 2D code to numerically integrate the system of equations (9) and 
(10), which describe the self-consistent evolution of the DW-ZF turbulence systems. We have 
chosen v in /uj ci = 0.1, va/ud = 0.01, o = 0.1, and p s /L n = 0.01. Numerical descritization 
employs the spatial derivative in Fourier spectral space, while time is descritized using time- 



221 ] . Periodic boundary conditions are 



split integration algorithm, as prescribed in Ref. 
used along the x and y directions. A fixed time integration step is used. The conservation of 
energy 23] is used to check the numerical accuracy and validity of our numerical code during 
the nonlinear evolution of the small scale drift wave fluctuations and zonal flows. We also 
make sure that the initial fluctuations are isotropic and do not influence any anisotropic 
flow during the evolution. Anisotropic flows in the evolution can, however, be generated 
from a k y = mode that is excited as a result of the nonlinear interactions between the 
ZFs and small scale DW turbulence. The ZF and DW fields are initialized with a small 
amplitude and uniform isotropic random spectral distribution of Fourier modes in a 2D 
computational domain. These fields further evolve through Eqs. (9) and (10) under the 
influence of nonlinear interactions. Intrinsically, the set of Eqs. (9) and (10) possesses 
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FIG. 1: Evolution of mode structures in our coupled DW-ZF turbulence model from an initial 
random distribution. In the presence of self-interaction terms, zonal flows are enhanced and quench 
the DW turbulence more efficiently. Numerical resolution is 256 2 , box size is 2ir x 2ir. 



parametrically unstable modes involving short scale drift waves and zonal flows. In the 
early phase of simulations, we obtain the growth of small scale DWs. We have carried out 
two characteristically distinct sets of simulations by switching on and off the self-interaction 
terms. This enables us to gain considerable insight into the physics of generation of zonal 
flows and associated transport level in the coupled DW-ZF turbulence systems. 

To gain insight into the characteristics nonlinear interactions in our coupled DW-ZF 
turbulence model, we examine the Hasegawa-Mima-Wakatani (HMW) model 20|, |21[ that 
describes the electrostatic drift waves in an inhomogeneous magnetoplasma. First, the ion 
polarization drift nonlinearity in the HMW model, viz. z x V0 ■ VV 2 0, signifies the self- 
interaction Reynolds stress that plays a critical role in the formation of the ZFs 12|. This 
nonlinearity is basically responsible for the generations of the ZFs. Secondly, since in the 
collisionless DW dynamics, the electron density perturbation follows a Boltzmann law due 
to the rapid thermalization of electrons along z, the nonlinearity z x • V0 comes from 
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FIG. 2: The self-consistent generation of zonal flows is shown. In the presence of the self-interaction 
nonlinearity, zonal flows are generated rapidly and their saturated level is also enhanced when 
compared with the evolution without the self-interaction nonlinearity. 

the cross-coupling of the ZF's Ej_ x z particle motion with the drift wave density fluctuation 
in our model. The role of this nonlinearity has traditionally been identified as a source of 
suppressing the intensity of the nonlinear flows in the DW turbulence 23|. Nevertheless, the 
presence of the linear inhomogenous background can modify the nonlinear mode couplings in 
a subtle manner. Our objective here is to understand the latter in the context of the coupled 
DW-ZF turbulence system. The initially isotropic and homogeneous spectral distribution 
associated with potential fluctuations, as described above, evolve dynamically following the 
set of Eqs. (9) and (10). 

The small amplitude initial drift wave fluctuations are subject to the modulational in- 
stability on account of their nonlinear coupling with ZFs. The parametrically unstable 
fluctuations grow rapidly during the early phase of the evolution. The instability eventually 
saturates via the nonlinear mode couplings in which the DW Reynolds stress, in concert 
with other nonlinearities in Eqs. (9) and (10), play a critical role. The mode couplings 
during the nonlinear phase of the evolution leads to the formation of non-symmetric zonal 
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FIG. 3: Evolution of the cross- field diffusivity in the presence (blue curve), as well as in the absence 
(red) of the self-interaction term. The cross-field diffusivity with and without zonal flows is shown 
by the dashed and solid curves, respectively. Clearly, the presence of the self-interaction term 
enhances zonal flows, which dramatically reduce the cross-field diffusivity. 

flow structures. This is shown in Fig. 1. Our simulations exhibit that the self-interaction 
terms not only suppress the modulational instability on a rapid timescales, but they also 
regulate the generation of the ZFs (see, Figs 1 b, c, e, f). The extent and amplitude of 
the ZFs in Figs. 1(b) and (e) [with the self- interaction] are larger than that of (c) and 
(f) [without the self-interaction]. The final (i.e. steady-state) structures, nonetheless, show 
the formation of a predominantly dipolar vortex in the DW fluctuations, while the ZFs are 
dominated by a large-scale monopolar-vortex motions. It is noteworthy that the absence of 
the self-interaction contaminates the flows with more small scale structures (See Figs. 1 c 
and f). 

We next investigate the quantitative evolution of the ZFs, which is depicted in Fig. 
2. The spectral transfer of energy in the ZF is estimated from X)fc \4>{kx,k y = 0)| 2 . The 
latter describes a pile up of energy in the k y = mode that is summed up over the entire 
turbulent spectrum. We find from our simulations that the presence of the self-interaction 
nonlinearity rapidly suppresses the linear phase of the modulational instability of the DWs. 
Hence, the linear phase during the evolution terminates rapidly when compared with the 
no-self-interaction case. Alternatively, the modulational growth rate is enhanced, and it 
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saturates on a rapid timescales. Consequently, the presence of the self-interactions gives rise 
to an enhanced level of ZFs, as shown in Fig. 2. 

A direct consequence of the enhanced ZFs is to markedly suppress the cross-field turbulent 
transport, because the sheared flows (in the poloidal direction) associated with the ZFs tear 
apart the DW fluctuations/eddies, and thereby keeping their amplitudes low. We have 
computed the cross-field diffusivity in our simulations by using the ion fluxes involving 
the Boltzmann electron density perturbation and the linear ion polarization drift velocity 
associated with the nonthermal DWs. The cross-field ion diffusivity reads 

where D B = cT e /eB is the Bohm diffusion coefficient and (fik is the spectral potential 
distribution of the DWs. Note that the ion thermal diffusivity, as defined above, is subdued 
in the presence of the ZFs, since the latter eventually suppress the DW turbulence so that 
the steady state cross-field transport level is reduced. Consistent with this scenario, we 
find, from our simulations, that the onset of the ZFs quenches the cross-field turbulent 
transport, as shown by the dashed-curve in Fig. 3. Furthermore, due to the vanishing 
poloidal wavenumber of the ZFs, the sheared flows do not cause any cross-field turbulent 
transport in magnetized plasmas. 

In summary, the most notable point that emerges from our simulations of the coupled 
DW-ZF turbulence system is the importance and significance of the self-interaction nonlin- 
earity in modeling the low-frequency DW turbulence that is believed to be a critical source 
for heat and energy losses. A most realistic and accurate understanding of the latter is, 
therefore, essential for the building and the performance of the next generation controlled 
thermonuclear fusion reactors, such as the ITER. In the work, we have, for the first time, 
brought about the importance of self-interaction processes and their role with regard to the 
cross-field turbulent transport in high-temperature plasmas of thermonuclear fusion devices 
(e.g. tokamaks). For our purposes, we used a new set of nonlinear equations for the coupled 
DW-ZF turbulence system, which is a generalization of Ref. 7, by including self-interactions 
among DWs which drives finite amplitude ZFs. Numerical simulations of the newly obtained 
nonlinear equations reveal that the coupled DW-ZF turbulence system evolves in the form 
of short-scale drift dipolar vortices and a large scale monopolar zonal flow structure. The 
simultaneous presence of the dipolar and monopolar vortices is responsible for a subdued 
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cross-field turbulent transport in a magnetically confined fusion plasma. 
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